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Mechanical Design, Modelling and Control of a Novel Aerial 

Manipulator 

Alexandras Nikou, Georgios C. Gavridis and Kostas J. Kyriakopoulos 


Abstract —In this paper a novel aerial manipulation system 
is proposed. The mechanical structure of the system, the 
number of thrusters and their geometry will be derived from 
technical optimization problems. The aforementioned problems 
are defined by taking into consideration the desired actuation 
forces and torques applied to the end-effector of the system. 
The framework of the proposed system is designed in a CAD 
Package in order to evaluate the system parameter values. 
Following this, the kinematic and dynamic models are developed 
and an adaptive backstepping controller is designed aiming to 
control the exact position and orientation of the end-effector in 
the Cartesian space. Finally, the performance of the system is 
demonstrated through a simulation study, where a manipulation 
task scenario is investigated. 

I. Introduction 

Aerial manipulation is a new scientific field which has 
been gaining significant research attention and a wide variety 
of structures have been proposed in the last years. These 
manipulation systems possess several features which have 
lately brought them in the spotlight, with their objective 
mainly oriented towards performing effectively complex 
manipulating tasks in unstructured and dynamic environ¬ 
ments. Having them include active manipulation as a major 
functionality, would vastly broaden the applications of these 
systems, as they move from mere passive observation and 
sensing to interaction with the environment. Therefore, new 
scientifically applicable horizons will be introduced related 
to cooperative manipulation, surveillance, industrial inspec¬ 
tions, inspection and maintenance of aerial power lines, 
assisting people in rescue operations and constructing in 
inaccessible sites by repairing and assembling. Naturally, 
both designing and controlling aerial manipulators could be 
considered as nontrivial engineering challenges. 

The first theoretical and experimental results on aerial 
robots interacting with the environment were developed 
in [1], [2] using a ducted-fan prototype UAV within the 
framework of AIRobots project. The design of a quadrotor 
capable of applying force to a wall maintaining flight stability 
was performed in [3]. In [4] experimental results with a small 
helicopter with grasping capabilities were derived, along with 
the stability proofs while grasping. Several grippers that al¬ 
low quadrotors to grasp, pick up and transport payloads were 
introduced in [5]. An implementation of indoor gripping 
using a low-cost quadrotor has been introduced in [6]. The 
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authors in [7] addressed the problem of controlling multiple 
quadrotor robots that cooperatively grasp and transport a 
payload in three dimensions. Another significant work with 
cooperative quadrotors throwing and catching a ball with a 
net was performed in [8]. 

A dexterous holonomic hex-rotor platform equipped with 
a six DoF end-effector that can resist any applied wrench was 
proposed in [9]. A system for aerial manipulation, composed 
of a helicopter platform and a fully actuated seven DoF 
redundant robotic arm, has been introduced in [10]. Another 
hex-rotor manipulator that consists of three pair of propellers 
with a two-link manipulator aiming to trajectory tracking 
control was studied in [11]. 

More recently, significant experiments using commercial 
quadrotors equipped with external robotic arms have been 
conducted in [12]—[14], 

In this work, a completely novel aerial manipulator is 
introduced. This could be considered as a small autonomous 
aerial robot that interacts with the environment via an end- 
effector by applying desired forces and torques in a 6 DoF 
task space. The proposed system provides mechanical design 
flexibility achieved through technical optimization problems. 
The structural geometric distribution is the outcome of the 
aforementioned problems with the main goals being oriented 
towards low body volume, controllability of the system, 
avoidance of possible aerodynamic interactions and effi¬ 
ciency in performing desired manipulation tasks in dynamic 
environments. The system is fully integrated as it is not a 
commercial aerial robot equipped with an external robotic 
arm, as many of the aerial manipulators mentioned above. 
The optimal number of thrusters, their positions/orientations 
and the optimal position of the end-effector on the body 
structure are defined with respect to the modelling design 
limitations. Taking all the above into account, the remaining 
challenge is is to actually construct this novel aerial robot. 

The rest of paper is organized as follows. In Section [H| a 
functional description of the robot and the mechanical design 
analysis is discussed. A mathematical model that captures 
the proposed system dynamics and govern the behaviour of 
the system is derived in Section [III] Based on this highly 
nonlinear model, an adaptive backstepping control law is 
designed in Section IV In Section [V] simulation results are 
presented in order to study the performance of the system. 


Finally, the main conclusions are discussed in Section VI 


II. Mechanical Design 

The overall description of the Aerial Manipulator was 
based on the idea of designing an aerial robot composed 




of a set number n of similar thrusters and an end-effector, in 
order to interact with objects in the environment. The exact 
geometry of the structure will be the result of the analysis 
of this section. 

A. Principles of the Problem 

Initially, we define the Body-Fixed frame and the End- 
Effector frame as F B = {x B ,y B , z B }, F E = {x E ,yE,z E }- 
These frames are attached to the rigid body of the aerial 
manipulator as in Fig. [I] The vectors r,, r e £ R 3 denote 
the position of each thruster and the position of the end- 
effector respectively with reference to the Body-Fixed frame. 
The thruster orientations are given by the unit vectors F. t £ 
R 3 , i = 1 the thrust forces are defined as A i and 

the propulsion vectors are given by A, F,. At this point, 
it is assumed that the total system is considered to be a 
rigid body and, without loss of generality, the End-Effector 
frame and the Body frame have the same orientation. Thus, 
the actuation force applied to the end-effector is F act \ B = 
F act \ E £ R 3 , where \ B ,\ E denote the expressions to the 
frames F B ,F B respectively. The corresponding actuation 
torque is obtained via the formula T act | = T act | + r e x f e 

where T act \ E = r x f e is the torque produced by the 
end-effector. The terms r, f e denote the displacement vector 
(length of the lever arm) and the vector force that tends to 
rotate a gripped object from the end-effector. 



Fig. 1: Aerial Manipulator Frame Configuration System 


B. Forces and Torques 

The forces transmitted essentially through the end-effector 
are written as 

n 

J2(X l F t ) + W = F act \ B (1) 

where W £ R 3 is the vector that corresponds to the total 
weight of the system. By separating the weight as w s = 
W — n w, where w £ R 3 is the weight of each thruster, 0 
is modified as 

F X + n w + w s = F act | B (2) 

where A = [Ai • • • X n } T £ R” and F=[F 1 --- F n ] T £ R 3xn . 


Similarly, the torque from each thruster is Tj = ry x 
(XjFj) = XiS(ri)Fi. The well-known skew-symmetric ma¬ 
trix S(-) £ R 3x3 is defined as a x b = 5(a) b for the 
cross-product x and any vectors a,b £ R 3 . The torque due 
to the weight is calculated as 

n n 

T w = rcxW = ^ ~2riXw+r s xw s = y ^ S(rj)w+S(r s )w s 

i =1 2=1 

where re is the centre of gravity of the system and r s 
is the centre of gravity of the system, when omitting the 
mass of each thruster (r// lhr ). The reaction torque of each 
thruster is n = p (A Ff) where p is a coefficient that 
represents the relationship between the thrust force and the 
reaction torque [15]. Therefore, by combining all torques the 
following equation holds 


E{a. s(n) F^ + p 53 (a 4 f.) 

2=1 2=1 
n 

+ '^2 s ( r i) w + S(r s ) w s = T act \ B (3) 


Using the matrices 


r=[n ••• r n \ £ R 3x ”,F = F act \ B - n w - w s (4) 
E(r,F) = [5(n) F 1 • • • S(r n ) F n \ £ R 3xn (5) 

in we get 

' F X = F 

( 6 ) 

o y vi) uj — o yr s ) w s 

i =1 

F 


E X = T act \ B - P F - 53 S(n) w - S(r s ) 


E(r,F) 


’, from 


By defining the matrix D(r,F) = 
the system |6]i it is implied that 

D(r, F) X = W R (7) 

where the augmented wrench vector W B £ R 6 is given by 

F 


W R = 


T ac t\ B ~ pF - 53 S{ri)w - S(r s )w s 


2=1 


( 8 ) 


C. Negative Thrust Forces 

It is clear that when solving ([7]), the vector that corre¬ 
sponds to the thrust force A can obtain any value in R 6 . 
However, the thrusters are optimally designed to produce 
thrust force towards a specific direction, which we set to 
correspond to the positive values of A,;. In order to alleviate 
the problem of negative A,, a conservative solution is adopted 
in this analysis, which is based on the idea of introducing 
one additional thruster. Thus, 0 is rewritten as 


53 Xi ti = Wr 

i —1 

where 

D{r, F) = [ti f 2 • • • t n \ and f = 


Ft 

S(n) Pi 


(9) 


( 10 ) 










for all i = 1,..., n. The vector 


n 

ta — ^ ^ ti 

2 = 1 


~Y.Fi 

2 = 1 


F a 

S{r a )F a 


that corresponds to the additional thruster is introduced. 
Using CLD’ the position vector r a and the direction F a 
of the new thruster should satisfy the equations F a = 

~Y.Fi, S(F a ) r a = -Y (sw) n\. 

If we assume that 0 results in some negative thrust 
forces, then the set oyv = {k : A k < 0, k = 1,6} denotes 
the indexes for every negative thrust force and crp = 
{1, 2, 3,4,5,6} — ctjv the corresponding set of positive thrust 
forces. Observing that Afc < 0 <£=> (—A k) >0, V k € oyv, 
0 can be separated in 


y^ Ai ti + yy Afc t k = wr 

2(E<T p k£(JN 

<=> ^2 ^ ** + yz (—Afc) (—tk) = Wr (12) 

i£crp k£(TN 

Now, from © the following can be exported 

n 

t a =—yy ti = — yy ti — yy tj 


iGo-jv 


(13) 


It is obvious that 

~ ^2 tj = ~tk - ^2 tj’ V k G ctjv (14) 

jGo-jv jGcTiv 

jAk 

Combining ( fl3| , ( fl4| we obtain 

- 4 = fq + yy u + yy fj, v k e o-jv d5) 

iGcr_p JGitjv 

j¥=k 

By substituting ( fl5] > into © we have 

yy Ai u + y^ (—A/c) f a + yy ^ + yy tj = wr 

i£(Tp ZcGcTiv L i^L<Jp j€&N J 

j¥=k 

Defining 


A = (—Afc) > 0, Ufc = yy (—Aj) > 0 (16) 

fcGo-iv j Gctjv 

j^k 

and rearranging the terms we result in 

y ( (Ai + A )ti + y ( E k t k + A t a = Wr (17) 

i£crp k£(TN 

From © the thruster redistribution among all thrusters after 
adding the new thruster is provided. It has been proven 
that the issue of negative thrust forces can be alleviated 
with adding one extra thmster. This equation can be better 
analysed in Fig.[2]in which the thrust redistribution algorithm 
is depicted. The variables A 7 , A 7 , X' k denote the initial thrust 
forces and the other variables the thrust force after the 


redistribution, plus the additional thrust force A a . Thus, the 
six thrust forces (not necessary all positive) are equivalent 
to seven thrust forces, all positive with redistributed thrust 
forces as in ©. By using the additional thruster, 0, 0 
are reformed into 


F = F act L - (n+ 1) w- 


W R = 


T a , 


F 

n+1 

— H F — yy S(ri) w - S(r 
2=1 


(18) 

(19) 



D. Aerodynamic Interaction 

At this point, the aerodynamic interaction between the 
operation thrusters is investigated. The aerodynamic effects 
produced by each thruster, are based on experiments that took 
place in Control Systems Lab NTUA on a 8 x 4.7SF APC 
propeller accompanied with the Neu Motor NEU 1902/2Y - 
2035 motor, which produces at 17550 rpm, a A max = 28 N 
thrust force. The surface, corresponding to every thruster, 
that approximates these effects is described by a third order 
equation in (SI). By expressing this equation in the thruster 
frame F T . = {i 7 , y\, z(}, i = 1,..., 6, a we get 

-0.06 < x'i < 0.91 (20) 

(y 7 ) 2 + (z'i? < [—l-l(a; 7 ) 3 + 1.56(:r 7 ) 2 - 0.3(x 7 ) + 0.11] 2 

Hence, the aerodynamic effects of the air flow throughout the 
rotor are extended from x = — 0.06m to x = 0.91?n. The x 
axis shows the length of the aerodynamic effect of the exit 
flow. In order to understand this, one should consider the 
rotor/blade to be positioned at x = 0. On the other hand, 
y axis shows the distance of the effect measured from the 
rotation axis of the blade, where at position (x = 0, y = 
0.102m) (SI) is the blade radius in approximation (because 
of the existence of an offset). 

An arbitrary point p = [x y z] T expressed in Fr and the 
corresponding p\ = [a; 7 y 7 z'] T expressed in Fp ., can be 
linked together by the equation p = {TR 1 ^ (fi, Fi)} p(, 

F'P , ±, , 

where TR F J (ay, F,) is the appropriate homogeneous frame 
transformation corresponding to the translation and orien¬ 
tation vectors (rj,P)). By combining the previous coordi¬ 
nates transformation equation with the constraints ©, a 
set of constraints that can be described in matrix form as 
G(ri,Fi,p ) < 0, is produced. The distance between two 
such volumes i, j can be defined and evaluated via the 
optimization problem (Pi) of 'Fable |TJ 
















E. Design Problem 

Given a particular structure defined by the matrices 
(r, F), for a set of required actuation forces and torques 
(Fact | Tact \ E ) it i s necessary to find the associate levels of 
the thrust forces A,. Since Wr G R 6 , in order for ([TJ) to have 
a solution for A G M™, the conditions {rank(D) = 6, n > 6} 
are required. The rank condition is adequate from a strict 
mathematical perspective but from a practical point of view, 
as 0 leads to the thrust forces values A G R n , the sought 
solutions should not be very sensitive to small deviations. 
This is partially achieved by using the condition number 
k(D) = (T max (Z))/<7 min (-D) where a(D) = ^eig (D T D) 
are the singular values of the matrix D, eig(-) denotes 
the eigenvalues of a matrix and cr max (-D), a m i n (D) are the 
maximum and minimum singular values of the matrix D 
respectively. Thus, a low condition number n(D) > 1 is 
required [16]. Although the condition number is bounded 
to take feasible values (not equal to zero/infinity) when 
a(D) —» 0, the matrix D(r,F ) might be ill-conditioned i.e. 
det(D(r, F)) —► 0. Thus, cr(D) > ei > 0. Furthermore, to 
avoid the fan interaction an other constraint is introduced 
as dij(ri,Fi,rj,Fj) > e 2 > 0, \/i,j = 1,2,... ,n, a. Note 
that, similarly to (Pi), the position r e should be introduced 
to the design problem as the intersection avoidance between 
a sphere (with radius R e ) that encloses the end-effector, 
and the thrusters. This sphere, when expressed in the End- 
Effector frame F E , is given by ( x' e ) 2 + ( y ' e ) 2 + ( z ' e ) 2 < R 2 e . 
Therefore, the constraint associated with the end-effector is 
d ei (r e , ri) > R e > 0, Vi = 1, 2,..., n, a. An optimization 
is also required to minimize the volume of the system, by 
using the norm J(r) = ||r'|| 2 ■ Taking all the above into 
consideration, the design problem is essentially recast to the 
optimization problem (P 2 ) from Table [I] The optimization 
parameters are chosen as K = 5, = 10~ 3 , e 2 = 

lCT 2 m, R e = 10 -2 to. 


(Pi) 

dij(ri,Fi,rj,Fj) = min | \p t - Pj\\ 

Pi iPj 

s.t. G(ri,Fi,pi)< 0 

G(ri, Fi,pj) < 0 

(P2) 

min J(r) 
r,r e ,F 

s.t. cr(D) > ei 

dij > € 2 , = 

d e i > R e , V i = 1,2 ,... ,n, a 

n 

i =1 

n 

S(F a ) r a = -^S(Fi) n 

1 < k{D) < R- 1 


TABLE I 

Optimization Problems 


F. Solving the Optimization Problem 

It should be noted that when solving the optimization 
problem (P 2 ), each time the inner problem (Pi) should be 
solved. There are 45 decision variables of the optimization 


problem, which correspond to the seven position vectors 
( ri ) of the thrusters, the position vector (r e ) of the end- 
effector and the direction vectors (P;) of the seven thrusters. 
This issue, entails the necessity of solving 28 optimization 
problems for each evaluation attempt of the outer problem 

0P 2 ). 

The inner problem, that refers to the avoidance of the fan 
interaction, is smooth but in terms of the outer problem 
(P 2 ) is nonsmooth and nonlinear. The objective function 
and the constraints of the problem (PI) are continuous 
and this problem, according to the inputs, has one and 
only one global minimum. Using the appropriate rotation 
and transformation matrices, the (PI) was solved by the 
active-set strategy [17],[18]. On the other hand, the design 
problem (P 2 ) has nonsmooth, discontinuous and nonlinear 
inequality constraints, but smooth objective function. Con¬ 
sequently, a non-gradient-based methodology that searches 
disjoint feasible regions, is utilized. For the pre-search of the 
design space, a Latin Hypercube (LHS) [19] was chosen, in 
order to ensure that the points are distributed throughout the 
search space. The Latin Hypercube sampling is known to 
provide better coverage than the simple random sampling 
[20]. Following this, a Generalized Pattern Search (GPS) 
direct search algorithm [21],[22] was used. 

The thrust force (A) and the momentum (Q) can be 
calculated from [23],[24] as 


Q = np Cq R 5 n 2 

A = 7 xpC x R 4 n 2 




Q=^R A 


( 21 ) 


where the term R corresponds to the coefficient //, R 
is the radius of the rotor and p. f l denote the air density 
and the rotational speed of the rotor respectively. Applying 
a combination of the Blade Element Theory [24] and the 
Momentum Theory [15], using the modified versions pro¬ 
posed in [23] and invoking the experimental results extracted 
by our lab on the APC propeller, it was calculated that 
C x = 0.008, Cq = 0.0095, pt = 0.1473, R = 0.124to. By 
solving the optimization problems, with the results depicted 
in Table |Tl] the matrix l)(r. F) is full rank and using (|7j, 
the thrust forces can be calculated as A = D~ 1 Wr. All the 
constraints were satisfied and a low volume body structure 
with condition number n(D) = 3.36 resulted. The wrench 
vector Wr can be determined by substituting the desired 
actuation forces/torques (F act | , T act | ) in m The maxi¬ 
mum thrust force and torque which can be applied from the 
system are A max = 28 N and 3Nm respectively. The values 
of the components, proposed for the Aerial Manipulator, 
are the following: the motor and the propeller (0.12 kg), 
the frame (0.66 kg), the battery (0.25 kg) and the electronic 
components (0.15 kg). The total mass of the proposed system 
is to = 1.90 kg. Ultimately, the production of a carefully 
studied framework (Fig. [3]) was achieved by using the 3D 
CAD Package SolidWorks. Using this Package, the system 
parameter values of the Table [Tl] have been evaluated. 







Fig. 3: Aerial Manipulator 3D Caption of the Framework 


Param. 

Description 

Value 

Units 

m 

Total Mass 

1.90 

kg 

m thr 

Thruster Mass 

0.12 

kg 

Ig 

Moment of 
Inertia Tensor 

0.3488 0.0683 - 0.0457 
0.0683 0.1588 0.0144 
- 0.0457 0.0144 0.4081 

kg m 2 

TG 

Centre of 
Gravity 
Position 

[0.0737 0.0083 - 0.0781] T 

m 

r e 

End-Effector 

Position 

[-0.23 0.015 0.23] T 

m 

r s 

Centre of Grav. 
from {3} 

[0.1267 -0.0052 -0.1900] T 

m 

J(r) 

Total Structure 
Volume 

1.80018 

m 3 

a 

Gravitational 

Acceleration 

9.81 

m/s 2 

n 

Thruster 

Positions 

n = [0.43 -0.15 — 0.44] T 
r 2 = [0.08 -0.22 -0.14] r 
r 3 = [0.1 - 0.9- 0.2] T 

r 4 = [-0.34 0.25 0.006] T 

r 5 = [0.184 0.359 -0.254] r 
rg = [-0.22 -0.44 -0.04] r 
r 7 = [0.51 0.79 - 0.06] r 

m 

Fi 

Thruster 

Orientations 

Fi = [0.08 0.39 0.92] r 

F 2 = [-0.33 - 0.90 0.29] r 
F 3 = [0.13 -0.87 — 0.48] r 

F a = [0.56 0.08 0.82] r 

# 5 — [0.83 0.11 — 0.55] T 

F 6 = [-0.66 0.57 - 0.49] r 
F 7 = [-0.59 0.62 - 0.51] T 



TABLE II 

Aerial Manipulator Parameters 


III. Mathematical Model of the Aerial 
Manipulator 

In this section, the kinematic and dynamic equations of 
motion in case there are no interaction forces and torques 
from the environment applied to the end-effector are pre¬ 
sented. 


and it should be noted that the Body-Fixed frame’s origin 
does not coincide with the centre of gravity G. The position 
of Fb relative to Fa can be represented by p = [x y z] T £ 
R 3 and the corresponding orientation by the rotation angles 
0 = [<fi 6 i/j] t £ R 3 . The translational and rotational 
kinematic equations of the moving rigid body are given (see 
[25]) in matrix form by 


V 


" MO) 

0(3x3) 

V 

© 


0(3x3) 

MO)_ 

UJ 


where 0( 3 X 3 ) is the 3x3 zero matrix, v = [v x v y v z ] T £ 
R 3 , u> = [u> x ujy u} z ] T £ R 3 denote the translational velocity 
and the angular velocity of Fb relative to F A respectively, 
both expressed in the Body-Fixed frame. The transformation 
matrices J t (0), J r (0) £ R 3x3 are given by 


MO) 


cgc^ 

S(j)SQCijj SqpC(p 

SQCfjyCip + SfiStip 

SxfjCg 

SjySoS^ + CfffCip 

SQS-iJjC^ SfjyCqp 

—sg 

S(f>C0 

C(j) CQ 


(23) 


MO) 


1 3^8 c<f>tg 

0 c<j) 

0 S(f>/CQ C(j) j CQ 


(24) 


The position of the end-effector with respect to F,\ is p e = 
[x e y e z e ] T = p+ MO) T e £ R 3 . Its derivative is obtained as 
p e = J t (0) v — Jt(0) S(r e ) w, using the formula J t (0) = 
Jt(Q)s{ui) from [26]. The Body-Fixed and the End-Effector 
frame have the same orientation with reference to Fa, as 
mentioned in Section [II] hence 0 e = 0. By combining the 
last results the following kinematic equation holds 


Pe 


'M&) 

-MO) S(r e ) 


V 

0. 


0(3x3) 

MO) 


UJ 


J(G) 


where 0 e is the orientation of Fe relative to F A . The 
Jacobian matrix of the system J(£ e ) £ R 6x6 relates in 
a straightforward way the linear velocity p e and the rate 
of change in the rotational angles 0 of the end-effector 
expressed in Fa, with the Body-Fixed velocities v, u>. 


B. Dynamic Model 

The dynamic equations can be conveniently written with 
respect to the Body-Fixed frame by using the Newton-Euler 
formalism (the main concept is discussed extensively in [27], 
[28]), as 

F 
T 

where 


M 


C(v) 



771/3 —mSira) 
mS(r g) Ib 


is the inertia matrix. 


M > 0, M = 0 


(27) 


A. Kinematic Model 

Fig. □ shows the reference frames defined to derive the 
kinematic and dynamic model of the proposed system. The 
Earth-Fixed inertial frame is defined as Fa = {i', 4 , y A , z A } 


C = 


mS(uj) —mS(uj)S(rG) 
mS(ra)S(u}) -S(I B u) 


, C = -C T (28) 


is the Coriolis-centripetal matrix, l A is the 3x3 identity 
matrix, m is the total mass of the system, Ib is the inertia 



















































tensor expressed in Fb and v = [v T uj t ] t £ R 6 is the 
vector of the Body-Fixed velocities. The inertia tensor can be 
written as Ib = I G — mS{r G )S(r G ) where Iq is the inertia 
tensor relative to the body’s centre of gravity. The vectors 
F\ ,T\ £ R 3 describe the forces and torques acting on 

the system expressed in the Body-Fixed frame and can be 
derived as 


PM 


F X — m g J[ (0) e 3 

PIbJ 


EX + pFX — mg S(r G ) J 4 T (0) e 3 


A + 


0(3x6) 

f,F 


X — m g 


h 

S(r G ) 


Jt ( 0) e 3 


propulsion 

forces/torques 


reaction 

torques 


gravitational 

forces/torques 


(29) 


where e 3 = [0 0 1] T . Combining 
respect to [v T u ; T ] T we get 


and solving with 


= —M~ 1 C(u) 


v 

UJ 

-1 




E + /j F 


— to g M 
= H(y) + G(p.) +N A 


S(r a ) ^ 63 


where the matrices are defined as 


H{v) = —M _i C(v) v,N = M 


-l 


F 

E + /i F 


(30) 

(31) 


> 0 (32) 


G(£ e ) = -m g M~ 


Jt (©) e 3 


S{r G ) 

B(£ e , v) = H{v) + G(p e ), B : R 6 x R 6 -) R 6 

IV. Nonlinear Control of the Aerial 
Manipulator 


(33) 

(34) 


matrix N is full rank with low condition number which 
constitutes a vital result of the control oriented optimization 
from Section HI1 

The system ( |35j ) is highly nonlinear, cascaded and fully 
actuated in the well-known strict feedback form, with vec¬ 
tor relative degree 2. For such systems, the backstepping 
controller design has proven to be successful [29], [30]. 
Due to the fact that the system is in the presence of the 
uncertainties 6\ and the disturbances d(£ e ,v,t), a robust 
adaptive controller will be designed in order to tackle them. 
The aim is to study if the proposed system with the resulting 
geometry from the optimization problems (Pi), (P 2 ), the 
system specifications from Table III and the aforementioned 
uncertainties/disturbances from ([35 1, is capable to perform 
specific trajectory tasks efficiently. In order to design the 
controller of the system ( |35j ), the following assumptions are 
required: 

Assumption 1: The states of the system p e , v are available for 
measurement Vi > 0 for the following control development. 
Assumption 2: The desired trajectories pd es are known and 
bounded functions of time (pj es £ C^) with known and 
bounded derivatives (pdes ) Pies £ Poo)- Assumption 3: The 
disturbance d(£ e ,v,t) = [di(£ e ,v,t) ■■■ de{S, e ,u,t)] T is 
unknown but bounded with |dj(£ e) v, t)\ < A * where A,; 
are unknown positive constants for all i = 1,..., 6 and 
t > 0. Assumption 4: It is assumed for all i > 0 that 
— \< 0(t) < This ensures that the Jacobian matrix is 
nonsingular since det(J(£ e )) = 1 /eg. This assumption is 
likewise utilized in [26], [31]. 

• Step 1: To begin with the backstepping controller design, 
the position-orientation error of the end-effector is defined as 
Z\ = £ e — pics £ R 6 - By differentiating it and using ( |25[ i we 
get 

Zl = J(£e) V - £des (36) 


A manipulation task is usually given in terms of the 
desired position and orientation of the end-effector. The 
objective of this section is to design a controller for the aerial 
manipulator ensuring that the position p e (t ) and the orien¬ 
tation 0(f) of the end-effector track the desired Cartesian 
trajectory pies(f) = [Pe,des(*) ^des^F e R6 asymptotically 
while all the closed loop signals remain bounded for all 
t > 0. Firstly, by using formulas 43T} the aerial 

manipulator model, including the kinematics and dynamics, 
can be written as 

ie — J(te) V 

y ’ (35) 

V = P(£e, v) + N 6\ A + d(£e, V , f) 

where d : R 6 x R 6 x R + —X R 6 represents the unmodelled 
nonlinear dynamics and the environmental disturbances. The 
unknown matrix 9\ = diag{0*,..., 9£} £ R 6x6 with 
9* £ [0 m i n , 0 max ] = [0.1,1], is introduced to model the 
control actuation failures and the modeling errors among 
the thrusters of the system, e.g. if 9* = 0.8 then the i —th 
actuator has 20 % controller effectiveness reduction. The 
control inputs of the system are the six independent thrust 
forces A,(f), i = 1,...,6 as mentioned in Section [S] The 


We view v as a control variable and we define a virtual 
control law z^es € R 6 for ( |36| . The error signal representing 
the difference between the virtual and the actual controls 
is defined as 22 = v — z/d es £ R 6 - Thus, in terms of 
the new state variable, ( |36| can be rewritten as i,\ = 
>A(£e) Z 2 + J(£ e ) z^des ~ £des- Consider now the positive 
definite and radially unbounded quadratic Lyapunov function 
V\{z\) = A|| zi || 2 = \z[zi- By differentiating it with respect 
to time yields 

Vj = z[il = zl j J(£e) I/des - pies} + ^J(Pe )^2 ( 37 ) 

The stabilization of Z\ can be obtained by designing an 
appropriate virtual control law 

I'des = J-\£,e) {Pies-i^i} ( 38 ) 

where the matrix Ki £ R 6x6 , Ki = Kl > 0 represents the 
first controller gain to be designed. Hence, the time derivative 
of Vi becomes V\ = — zlK\Z\+zl J{^ e )z 2 - The first term of 
on the right-hand of this equation is negative and the second 
term will be canceled in the next step. 

• Step 2: For the second step we define the matrices of the 
parameter estimation errors as A = [Ai • • ■ Ag] T = [(Ai — 





























At) • • -.(Ae-A 6 )]t and 9 X = diag {(ft-fif), , (0 6 -^)} 

where A,;, 0, are the estimations of the unknown parameters 
A i,8* respectively. The time derivative of the error z 2 is 
i 2 = B(£ e , v) + N 9\ A + d(£ e , v, t) — t>des- The Lyapunov 
function candidate in this step is chosen as 

Vi{zi,Z2, A, h) = + z 2 +- A r r A 1 A+-tr(0^r e 1 8 x ) 

where 1’^ =Tg > 0.Fa = fj\ > 0 are diagonal adaptation 
gain matrices and tr(-) denotes the matrix trace. The time 
derivative of V-iizi, z 2 , A, 9 X ) is obtained as 

V 2 = -zlK lZl + zj { J T {£,e)Z\ + B(£e, v) + N 6\ A — f'des} 

+ a£d(&, V,t) + A r r^ 1 A + tr^r- 1 ^) (39) 

Using the -z\Ki Zl < — A mi „(^ 1 )||-^l|| 2 , z\ d(£ e ,v,t) < 
zJsgn(^ 2 )A and adding and subtracting the terms 
z^N9 x A, zJsgn( 2 2 )A in ([39} the following inequality 
holds 


V 2 < -A m in(A'i)||zi|| 2 + Z%{j T (€e)zi + B(£ e ,v) - i>des 
+ sgn(z 2 ) A + N 9 X \} - Z2N9 x A 
- 4sgn(* 2 ) A +A r rx 1 A + tr(^r-^ A ) (40) 

where \ m i n (Ki) denotes the minimum eigenvalue of matrix 
Ki, sgn (z 2 ) = diag{sgn(z 2 ,i),...,sgn(z 2i6 )} and sgn(-) 
denotes the sign function. Rearranging the terms and using 
the property a T b = tr(6 a T ), Va, b £ R n we get 


V 2 < —A m in(-^ 1 )||^l|| 2 + Z 2 {J T {£,e)Z\ + B(£ et v) — 2>des 

+ sgn( 2 2 ) A + TV 9 X A} + A r {r^ 1 A - sgn(z 2 )z 2 } 

+ tr{6 T x (T0 1 § x -N T z 2 X T )} (41) 

Given the form of V 2 from ( pT) the adaptive control law 
and the corresponding parameter estimation update laws for 
the nonlinear system ( [35] ) to be designed, are 

A(£e, v, A, 9 X ) = (§ X )- 1 (V- 1 {f'des - B(£ e ,v) - J T {^)z X 
- sgn(^ 2 ) A - K 2 z 2 } (42) 

A = r A {sgn(z 2 )z 2 - ° A} (43) 

k =T e Pto)(8 x ,N t z 2 X t ) (44) 

where K 2 = K T Z > 0 is the second controller gain matrix, 
a is a strictly positive gain (cr-modification rule [32]) and 
the projection operator Proj(-,-) is the same as the one in 
[33] with the parameter 6 to be designed. By substituting 
( |4 2 [ ) , ( |43) , ( |44| ) into ( 141) and using the property A T A = 

-||A|| 2 + -||A|j 2 — —1|A|| 2 the following inequality holds 


V 2 < — A m i n (ifi)||2l|| 2 — Amin (K 2 )\\z 2 I 


tr \e\ 


Proj(0 A! y)-yj)-J||A|| 2 -A||A|| 2 . 


I A|| : 


<0, y = N T z 2 \ T <-§||Ap+f||Ap 

The projection operator invoked from [33] contributes to the 
negative semi-negativeness of the Lyapunov function since 


by definition tr 


{^a Pr °j(^A ,y)-y } 


< 0, Vy. Moreover, 


it guarantees that if (9^(0) £ [^min^max] is chosen, then 
§i(t) £ [^min — 8, (9 max + <5], Vi = 1, 6, Vi > 0 for suitable 
S > 0. The last result protects the term ((9a) -1 in ( |42) from 
singularity. By defining w = |||A|| 2 > 0 we result in 


< -A m i„(iTi)||0i|| 2 -A min (iT 2 )||z2|| 2 -|||A|| 2 +«) 


from which it follows that both errors z\. z 2 and the param¬ 
eter estimation A are uniformly ultimately bounded with re¬ 
spect to the sets Ui = jzi £ R 6 : ||zi|| < ^/u)/A m i n (A'i) j, 

^2 = jz 2 e R 6 : ||z 2 || < yj fc/A m in(A 2 )^ and n A = 

| A £ R 6 : ||A|| < a/ 2 w/crj. Invoking that z\,z 2 are 
bounded and £des>Uies £ C^ then £ e ,i/ £ £<*,. Since 
A, A, 9 Xl $a, J^des are bounded then A, 9 X: A £ Uoq. Based 
on the above, it is proven that all closed loop signals remain 
bounded. 

One important issue associated with the controller de¬ 
sign is the analytical form of the time derivative of 
Vdes, which can be obtained from ( [38] ) as z>d es = 
J -1 (Ze)|/de S -i(£e) t'des — K\ z\ >, and the time deriva¬ 


tive of J(/ e ), which can be calculated by using the J r (0) = 


; 9 Jr „ N 

TT </>+ 0 , J{£e) = 


dcj) 
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J t (0) —J t (Q)S(r e ) 

Q(3x3) Tr(0) 



* [sac] 


Fig. 4: Position and Orientation Errors 


Required Thrust Forces 
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Fig. 5: Required Thrust Forces Using the Algorithm ([17} 


V. Simulation Results 

In this section, the results of a numerical simulation sce¬ 
nario are presented in order to demonstrate the performance 
of the proposed system. The dynamic model in ([35} is uti¬ 
lized with system parameters which are depicted in Table [TT] 




















































15 % controller effectiveness reduction is chosen with 9* = 
0.85 diagjl, 1, 1, 1,1,1}. The end-effector is forced to track 
the trajectory p e ,d es (t ) = [cos(0.5f) sin(0.5i) 1.5 + 0.3f] r 
with regulated orientation at @d es = [f f — f ] with 
reference to the Earth-Fixed frame. The initial conditions 
of the system are set to p e (0) = r e , p( 0) = 0(0) = 
v(0) = w(0) = 0(3 X i). The disturbance is set as d(t) = 
[0.5 0.4sin(2f) 0.4cos(f) 0.5 0.5cos(0.8f) 0.6sin(f)] T rep¬ 
resenting the unmodelled forces/torques and external distur¬ 
bances. The initial values of the parameters estimations are 
set to f?i(0) = ••■ = 0 6 (0) = 0.7 and Ai(0) = ••• = 
Af;(0) = 0. The controller gains are chosen as K\ = 
diag{l, 1,1,0.3,0.3,0.3}, if 2 = 8 diag{l, 1,1,1,1,1}. The 
adaptation gains are selected as a = 1.5, Ta = 
13 diag{l, 1,1,1,1,1}, T e = 0.1 diag{l, 1,1,1,1,1}. The 
parameter of the projection operator is set to 5 = 0.05. Fig. [4] 
shows the position and orientation tracking errors. The thrust 
forces are provided in Fig. [5] This paper is accompanied by a 
video demonstrating the simulation procedure of this Section. 
Due to space limitations, a video with an additional Scenario 
in better quality (HD) can be found at 

https://www.youtube.com/watch?v=DXnzu6XOrXs 


VI. Conclusions 

Aerial robots physically interacting with the environment 
could be very useful for many applications. In this paper, 
we have presented the mechanical design of a novel aerial 
manipulator which was the result of technical optimization 
problems. A mathematical model for the kinematics and dy¬ 
namics was derived in order to design an adaptive nonlinear 
controller to study the system while performing manipulation 
tasks. The simulation results illustrate the effectiveness of 
the proposed system and the controller to achieve tracking 
irrespectively of actuator failures, unmodelled dynamics and 
external disturbances. Future work mainly involves the con¬ 
struction of the aerial robot and the conduction of experimen¬ 
tal trials for the proposed framework with the actual system, 
in order to verify the theoretical results of this paper. 

References 

[1] L. Marconi, R. Naldi, and L. Gentili, “Modelling and Control of a 
Flying Robot Interacting with the Environment,” Automatica, vol. 47, 
no. 12, pp. 2571-2583, 2011. 

[2] A. Keemink, M. Fumagalli, S. Stramigioli, and R. Carloni, “Mechani¬ 
cal Design of a Manipulation System for Unmanned Aerial Vehicles,” 
IEEE International Conference on Robotics and Automation (ICRA), 
pp. 3147-3152, 2012. 

[3] A. Albers, S. Trautmann, T. Howard, T. Nguyen, M. Frietsch, and 
C. Sauter, “Semi-Autonomous Flying Robot for Physical Interaction 
with Environment,” IEEE Conference on Robotics Automation and 
Mechatronics (RAM), pp. 441^4-6, 2010. 

[4] P. Pounds, D. Bersak, and A. Dollar, “Grasping from the Air: Hov¬ 
ering Capture and Load Stability,” IEEE International Conference on 
Robotics and Automation (ICRA), pp. 2491-2498, 2011. 

[5] D. Mellinger, Q. Lindsey, M. Shomin, and V. Kumar, “Design, Mod¬ 
eling, Estimation and Control for Aerial Grasping and Manipulation,” 
IEEE/RSJ International Conference on Intelligent Robots and Systems 
(IROS), pp. 2668-2673, 2011. 

[6] V. Ghadiok, J. Goldin, and W. Ren, “Autonomous Indoor Aerial 
Gripping Using A Quadrotor,” IEEE/RSJ International Conference on 
Intelligent Robots and Systems (IROS), pp. 4645^-651, 2011. 


[7] D. Mellinger, M. Shomin, N. M. Nathan, and V. Kumar, “Coopera¬ 
tive Grasping and Transport Using Multiple Quadrotors,” Distributed 
Autonomous Robotic Systems, vol. 83, pp. 545-558, 2013. 

[8] R. Ritz, W. Muller, M. Hehn, and R. D’Andrea, “Cooperative 
Quadrocopter Ball Throwing and Catching,” IEEE/RSJ International 
Conference on Intelligent Robots and Systems (IROS), pp. 4972-4978, 
2012 . 

[9] G. Jiang and R. Voyles, “Hexrotor UAV Platform Enabling Dextrous 
Interaction with Structures-Flight Test,” 2013 IEEE International 
Symposium on Safety, Security, and Rescue Robotics (SSRR), pp. 1-6, 
2013. 

[10] F. Huber, K. Kondak, K. Krieger, D. Sommer, M. Schwarzbach, 
M. Laiacker, I. Kossyk, S. Parasel, S. Haddadin, and A. Albu-Schaffer, 
“First Analysis and Experiments in Aerial Manipulation Using Fully 
Actuated Redundant Robot Arm,” (IROS), 2013. 

[11] M. Kobilarov, “Nonlinear Trajectory Control of Multi-Body Aerial 
Manipulators,” Journal of Intelligent and Robotic Systems, vol. 73, 
no. 1-4, pp. 679-692, 2014. 

[12] S. Kim, S. Choi, and H. Kim, “Aerial Manipulation Using a Quadrotor 
with a Two DOF Robotic Aim,” IEEE/RSJ International Conference 
on Intelligent Robots and Systems (IROS), pp. 4990-4995, 2013. 

[13] M. Orsag, C. Korpela, and P. Oh, “Modeling and Control of MM- 
UAV: Mobile Manipulating Unmanned Aerial Vehicle,” Journal of 
Intelligent and Robotic Systems, vol. 69, no. 1-4, pp. 227-240, 2013. 

[14] A. Jimenez-Cano, J. Martin, G. Heredia, A. Ollero, and R. Cano, 
“Control of an Aerial Robot with Multi-Link Arm for Assembly 
Tasks,” (ICRA), pp. 4916^921, 2013. 

[15] G. Padfield, “Helicopter Flight Dynamics, The Theory and Application 
of Flying Qualities and Simulation Modelling” . AIAA education 
series, American Institute of Aeronautics and Astronautics, 1996. 

[16] L. N. Trefethen and D. Bau, “Numerical Linear Algebra”. SIAM, 
1997. 

[17] P. Gill, W. Murray, and M. Wright, “ Numerical Linear Algebra and 
Optimization”. Addison-Wesley Publishing Company, 1991. 

[18] K. Murty, “Linear Complementarity, Linear and Nonlinear Program¬ 
ming”. Sigma Series in Applied Mathematics, Berlin: Heldermann 
Verlag, 1988. 

[19] M. Stein, “Large Sample Properties of Simulations Using Latin 
Hypercube Sampling,” Technometrics, vol. 29, no. 2, pp. 143-151, 
1987. 

[20] M. McKay, R. Beckman, and W. Conover, “A Comparison of Three 
Methods for Selecting Values of Input Variables in the Analysis of 
Output from a Computer Code,” Technometrics, vol. 42, no. 1, pp. 55- 
61, 2000. 

[21] C. Audet and J. Dennis, “Analysis of Generalized Pattern Searches,” 
SIAM Journal on Optimization, vol. 13, no. 3, pp. 889-903, 2003. 

[22] C. Audet and J. Dennis, “A Pattern Search Filter Method for Nonlinear 
Programming without Derivatives,” SIAM Journal on Optimization, 
vol. 14, no. 4, pp. 980-1010, 2004. 

[23] V. Martinez, “Modelling of the Flight Dynamics of a Quadrotor 
Helicopter,” Master Thesis, Cranfield University, 2007. 

[24] R. Prouty, “Helicopter Performance, Stability and Control”. R.E. 
Krieger Publishing Company, 1995. 

[25] G. Antonelli, “Underwater Robots”. Springer Tracts in Advanced 
Robotics, Springer International Publishing, 2013. 

[26] T. Madani and A. Benallegue, “Backstepping Control for a Quadrotor 
Helicopter,” IEEE/RSJ International Conference on Intelligent Robots 
and Systems (IROS), pp. 3255-3260, 2006. 

[27] T. Fossen, “Guidance and Control of Ocean Vehicles”. John Wiley 
and Sons, 1994. 

[28] D. Bernstein, “Geometry, Kinematics, Statics and Dynamics”. Uni¬ 
versity of Michigan, 2013. 

[29] M. Ki'stic, I. Kanellakopoulos, and P. Kokotovic, “Nonlinear and 
Adaptive Control Design”. Wiley, 1995. 

[30] J. Zhou and C. Wen, “Adaptive Backstepping Control of Un¬ 
certain Systems: Nonsmooth Nonlinearities, Interactions Or Time- 
Variations”. Springer, 2008. 

[31] M. Huang, B. Xian, C. Diao, K. Yang, and Y. Feng, “Adaptive Track¬ 
ing Control of Underactuated Quadrotor Unmanned Aerial Vehicles 
via Backstepping,” American Control Conference (ACC), 2010. 

[32] E. Lavretsky and K. Wise, “Robust and Adaptive Control: with 
Aerospace Applications”. Springer, 2012. 

[33] H. Khalil, “Adaptive Output Feedback Control of Nonlinear Systems 
Represented by Input-Output Models,” IEEE Transactions on Auto¬ 
matic Control, vol. 41, pp. 177-188, Feb 1996. 


